Low energy modes and Debye behavior in a colloidal crystal 
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Abstract 

We study the vibrational spectrum and the low energy modes of a three dimensional colloidal 
crystal using confocal microscopy. This is done in a two-dimensional cut through a three- 
dimensional crystal. We find that the observed density of states is incompatible with the standard 
Debye form in either two or three dimensions. These results are confirmed by numerical simula- 
tions. We show that an effective theory for the projections of the modes onto the two-dimensional 
cut describes the experimental and simulation data in a satisfactory way. 

Keywords: Colloids, Phonons or vibrational states in low-dimensional structures and nanoscale 
materials. Lattice dynamics - Measurements. 



1. Introduction 

Colloidal particles dispersed in a solvent form fluid and solid phases under appropriate con- 
ditions . The ordered phase is interesting due to their long-range order combined with very 
soft mechanical properties making thermal fluctuations very important. For the colloidal systems 
considered here, the shear modulus is only a few Pa, whereas crystalline solids typically have 
moduli in the GPa rangefl Ql- Comparison of the properties of colloidal systems with harder 
molecular solids is a matter of ongoing research [5, 6, 7|. Earlier studies have largely focused 
on the thephonon dispersion behavior measured by means of video microscopy Jsl] or light scat- 
tering ['s'.'q']. However, to our knowledge there is no experiment that directly verifies the Debye 
scaling in the measured density of states of the vibrational modes. Even for molecular crystals, 
the this has proven rather hard to measure directly; mostly neutron scattering experiments have 
provided some but not much data that agree with the expected behavior flO'l. Perhaps for this 
reason the temperature dependence of the specific heat is usually taken as the hallmark for the 
Debye behavior of normal crystalline solids. The question we ask in this paper is how the spec- 
trum measured in a colloidal crystal compares with the predicted Debye behavior, DioS) ~ o)''^^ . 
To study the energy spectrum, we calculate the normal modes of a hard sphere colloidal crystal 
from the correlations in particle displacements. The colloidal particles used are small enough so 
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that they perform small 'vibrations' in response to the thermal excitations of the surroundings. 
We use confocal microscopy to observe and record these particle motions in a two dimensional 
field of view. The eigenvalues of the dynamical matrix and normal modes are then obtained from 
the spatial correlations ^ [THof the measured displacements. The present experimental results 
are compared and complemented with a simplified continuum theory as well as Monte Carlo 
simulation of hard sphere crystals. The main conclusion is that the spectrum of slice of a three 
dimensional system has a anomalous density of states which varies as D(a)) ~ o)^, which con- 
trasts the assumption made in [IZl that a 2d cut should exhibit the expecred 3d Debye behavior 
However our experimental results do agree with simulations presented below, and can in addition 
be explained by a simplified continuum theory. 

1.1. Experiments 

Hard sphere colloidal systems undergo phase changes with the volume fraction as control 
parameter There is no liquid-gas transition, but fluid-solid coexistence HI 01 is observed from 
the freezing transition point (pf - 0.494 until the melting point 0,„ - 0.545. A stable crystalline 
phase exists above 0.545 until the closed packed density at <p - 0.74. The colloids we use are 
charge stabilized PMMA particles with a diameter of a = 1 .3/im with a very small size poly- 
dispersity of about 2%. The particles are dyed with rhodamine and are suspended in a CHB 
(cyclohexyl bromide) / decalin mixture which closely matches both the density and the index of 
refraction of the particles. We further add organic salt TBAB (tetrabutylammoniumbromide) to 
screen any possible residual charges. 

The crystal was grown in a sample cell made of parallel plates with a confinement of approxi- 




Figure 1 : Snapshot of a two dimensional section in a field of view of 60fjm X 60//m with a hexagonal symmetry. 

mately ~ 1mm along the vertical direction. The volume fraction of the present colloidal crystal 
is about cf> x 0.57. Using confocal microscopy we acquire images at a speed of 25 frames per 
second of a two dimensional section of about 60 yum x 60 fim of the larger three dimensional 
crystal. The 2D slice was taken at a distance of 25 - 30 micron away from the coverslip, deep 
enough to avoid the effects of confinement. The entire crystal is polycrystalline, but we take our 
data from a region of the crystal that as far as we can see is perfectly crystalline and contains no 
defects. The particle positions are identified and tracked for a period of about 120 seconds using 
standard particle-tracking software; this results in a total of 3000 frames. Fig. [T]shows a typical 
snapshot of a two dimensional section of the measured crystal. 
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Figure 2: Left: Density of states as obtained from the correlations, eq.|4] Rigiit: The dispersion curves of transverse 
and longitudinal vibrations. The frequencies (vertical axis) near the zone boundary where dw/dk a coincide with the 
positions of the singularities (peaks) in the density of states on top. 



1.2. Density of states 

We obtain an ensemble of projected particle positions r,- - {xi,yi} from the above measure- 
ments. The displacement components from the mean position for the i'^ particle is given by, 

M/i(0 = r^i - <'"///>■ (1) 
Uv(i) = ryi - (ryi). 

where /i, v = x,y and "()" indicates an ensemble average. The k''' Fourier component of the 
above displacements is given by, 

N 

Uf,(k) = ^ Uf,(n)exp{-jlc.ri) (2) 

/■=! 



Uv(k) = ^ Uy(n)expi-jk.n), k = {k^,, ky) 



The summation above runs over all the = 1386 colloidal particles present in the system and 
the respective k values are chosen from the first Brillouin zone of the experimental lattice. Now, 
the potential energy of a harmonic crystal can be written as, 

U ^\/2j^{u;{k)K^y(k)Uym (3) 

k 

where K^^vik) fstlis the 2x2 dynamical matrix in Fourier space. From equipartition each of the 
above quadratic terms u*^{k)Kf,y(k)uy(k) /2 contains an energy of ksT/l. Therefore we can write, 

K^i(k)^{u;(k)uy(k))/kBT. (4) 
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In the present study, we take our two-dimensional data and analyse it according to the above 
scheme. The mode frequencies of the system are related to the eigenvalues A of the matrix 
^^1Was, 

cj(k) = ^Jl/Aik). (5) 

Now, there are of such 2x2 matrices (K~j) corresponding to the number N of k values. 
Separate diagonalization of each of these gives us 2N normal mode frequencies. Each of the 
branches of frequencies as a function of k then leads to the two dispersion curves shown in 
Fig. 12] The density of states D(a)) as obtained from these frequencies is shown in the same 
figure. A two peak structure is apparent in the density of states. Such peaks are familiar from 
the theory of lattice dynamics and often signal van Hove singularities |10] occurring due to 
the vanishing group velocities Vjjw at certain wave vectors. Indeed we see that the frequencies 
where doj/dk « (close the zone boundary) on the dispersion curves coincide with the ones 
around which peaks appear in the density of states. 

One important additional remark is that if we make the reasonable assumption that the system 
is in local equilibrium, the values of the squared-displacements are independent of the dynam- 
ics, which may be overdamped or purely ballistic, and even contain hydrodynamic interactions. 
What depends on the nature of the dynamics is the actual time-dependence of the displacements, 
but not their statistical distribution. Our new method [11, 14] directly analyzes the normal modes 
of the vibrations in the cages from the displacement correlations. If the system were harmonic, 
undamped and without hydrodynamic interactions, the DOS obtained from e.g Fourier trans- 
forming the velocity autocorrelation function and our method would coincide. However, in the 
former method all the spatial in- formation is lost. The strength of the present analysis using 
confocal measurements lies in directly visualizing the collective modes at low spatial frequen- 
cies and therefore gaining information about the nature of such modes. But the price to pay for 
this is the exclusion of knowledge of damping, anharmonicity and hydrodynamic interactions. 

1.3. Visualization of modes 

The next question we ask is what the low energy modes of the above system look like ? To 
be able to visualize these modes we compute the spatial correlation between the particles in real 
space via the covariance matrix. This matrix, which has been used recently in studying normal 
mode properties of colloidal LSi [Isl [I2I [HI 13 and granular systems Iil5i] . is defined as; 

Qttvii, j) = {u^,{i:))uy{ j)^ . (6) 

This is a 2N x 2N dimensional matrix for particles. Any eigenvector (v/) of the above ma- 
trix represents a normal mode at a single "frequency" aii = y/TJTi. In our colloidal system, it 
allows us to compute the normal modes, rather than supposing plane waves as in eq. (|4|. A few 
examples of the normal modes are shown in Fig. Q. The modes in the low-frequency part of 
the spectrum show a clearly coherent motion extending over large part of the field of view. For 
higher frequencies the modes appears to be rather random in nature. 
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Figure 3: Tomography of the normal modes shown in a field of view of 60/.i;ji X hO/mi. The pictures correspond to 
increasing frequencies going from a) to b) to c). The very low frequency ones show the extended plane wave nature of 
the modes. This coherent character is lost and modes appears to be random as higher frequencies are approached. 



1.4. Scaling of the low frequency density of states 

Having described the vibrational spectrum of the experimental crystal we come back to one 
of the central questions of this paper, which is what happens if one tries to study the matrix 
eq. Q using only a subset of the particles of the true three dimensional system, in particular 
if one reconstructs the dynamical matrix from the measurements from a single plane within the 
sample as we do. What is the nature of the effective interactions found in the analysis and 
what is the density of states found in the two-dimensional slice ? To answer this we look at 
how the measured cumulative density of states of the colloidal crystal compares to the expected 
integrated density of states NioS) - J D{(jj)d(jj ~ ii? \n d - 3 We find that the measured NioS) 
scales approximately as w"* rather than 0? . We fit the integrated density of states to a power law 
with adjustable exponent, a. Our data best fits to a value close to a = 3.9 + / - 0.19, far from 
the expected value a = 3 for a three dimensional elastic medium and even further from the value 
a - 2 expected in two dimensions. The quoted error is the one-sigma deviation of the Gaussian 
distribution of the scatter of the data around the fit. A standard error calculation gives that the 
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Figure 4: A power-law fit to the low frequency part of the cumulative density of states N{oj). The measured exponent 
(N(cj) ~ cj"), is a » 3.93 ±0.19. 



99 percent probability interval is 3.93 + / - 0.416, excluding the power 3. Furthermore, we have 
done a few different measurement series on the same sample, but at different locations within the 
sample and size of the region probed. Each of these measurements was analyzed separately to 
ensure the robustness and consistency of the results. Taking data from different regions in the 
sample, the exponent was found to always lie in the interval [3.7-4.1], with a mean value of 3.9. 
This again confirms the value of the exponent and its error quoted above. 

The above scaling of the low frequency cumulative DOS is also checked with the spectrum 
from the spatial covariance matrix that contains the non-diagonal correlationsfsec lZTt . This 
yields an exponent of about 3.6 -H /-0.3 which is reasonably close to the scaling (3.93 -H / -0.19). 
Thus, within the experimental error, the scaling of the DOS from both the full diagonalization 
and K^y are the same, and different from the usual Debye scaling. This new result will be 
further discussed and explained below, using Monte Carlo simulations of hard sphere crystals 
and arguments from elasticity theory. 

2. Theory and Simulation 

We check our experimental findings using Monte-Carlo simulation of a hard sphere crystal. 
We use a face-centered cubic crystal made of 864 particles, with periodic boundary conditions. 
We thermally equilibrate the system over a sufficient number (10^) of MC steps. We compute the 
covariance matrix (Eq. |6) to find the eigenmodes of the system along runs of 10^ MC sweeps, out 
of which we use 5x10"^ snapshots for time averaging. We obtain the matrix for both the whole 
three dimensional system and two dimensional (111) planes, in order to mimic the experimental 
conditions (see also [141 for a comparison between the 3d simulations and a 2d cut). As the 
number of particles (and thus modes) in a two dimensional plane is quite small, we average over 
several planes. The density of states we obtain is shown in Fig. |5] For the low frequency part 
of the cumulative DOS of the two dimensional slice we again find a cumulative density of states 
which seems to converge to N{cij) ~ of" when we increase the system size, which is compatible 
with the experimental data. Data with - 32000 can be fitted by exponents between 3.7 and 
4.1." 
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Figure 5: Left: Density of states for a two dimensional section of a 3D liard spliere crystal (with N = 4000) from Monte- 
Carlo simulation. Right: Low frequency part of the cumulative DOS for two dimensional slices of the thi'ee dimensional 
system for three system sizes, with logarithmic scale, showing a powerlaw behavior compatible with the Debye behavior 
N(ijj) ~ Lij* . Note that the range over which the scaling holds increases with the system size, indicating a convergence to 
the expected result. Frequencies are in unit of the inverse particle diameter. 



We now turn to a simple scalar theory which gives an indication as to how Debye theory must 
be modified when one observes two dimensional cuts of a three dimensional sample. 

The elastic energy of a three dimensional crystal can be written in terms of a symmetric 
shear tensor Ujj and three independent elastic constants. At large length scales correlations in 
displacement fluctuations decay as l/r, but are also characterized by a complicated tensorial 
structure coming from the cubic anisotropy of the crystal. 

The effect of the large distance decay can be found in a much simpler theory based on a 
scalar field u rather than the vector m,-. This scalar can be thought of as being, for instance, the 
amplitude of longitudinal fluctuation which couple to density fluctuations. The advantage of such 
a description is an enormous simplification in the tensorial algebra and a simple closed form for 
the projected correlation function. It is possible to perform a detailed tensorial calculation, which 
we will publish in the future and which yields very similar results. 

We therefore consider fluctuations of a scalar quantity u with an energy which is of the form 

f/ = ^ J (V«)'«f'r (7) 

where A is an elastic modulus. In the generalization to elastic fluctuations one would consider 
an energy based on the symmetrized strain tensor In Fourier space the energy has the form 

= (8) 

k 

we notice the usual scaling of the elastic energy in k^. 

In an underdamped system with kinetic energy piP'/l this gives rise to the dispersion relation 
aP' - c^k^. One thus expects a density of states 

dN ~k^Ak~ (J-Au (9) 

It is this scaling of the density of states in that is known from the theory of Debye. 
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In a three dimensional sample thermal fluctuations excite the system, so that equipartition 
and eq. (O implies 

<I«.P> = ^ (10) 

This gives a decay of correlations in real space which is given by the inverse Fourier transform 
of 1 Ik^. We can find the result immediately by reference to electrostatics: 

<mO>0-)> = — ^ ■ (11) 

AnA\ri - rj\ 

a Coulomb like decay of correlations. 

Now take a two dimensional slice of the system. Within this slice the correlations are still 
decaying as 1 /r. We wish to describe what we see, however, in terms of a purely two-dimensional 
theory, so we perform a two dimensional Fourier transform to find the efifective stiffness. Thus 
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where we use the subscript 2 to indicate that we are working with the two-dimensional projected 
objects. The result is rather interesting: rather than correlations in three dimensional being de- 
scribed by a decay in Ijk^ we find a slower decay: 1 /|fc2l- 

Now that we have the scaling form of the correlations we can work backwards and deduce 
the effective elastic theory in two dimensions 

f/2 = ^y'2feiKi2 (13) 

Thus the elastic behavior in real space corresponds to fractional derivatives of the field u leading 
to long-ranged effective interactions in the projected system. 
We now calculate the "propagative" eigenvalue by defining 

u\ = lAM (14) 

which is the analogy of uP- - c^k^ that we use in three dimensions. We note that the disper- 
sion law is very different from that of usual elastic problems. The density of states of this two 
dimensional matrix are just 

dA?2 ~ k2dk2 ~ -^dw2 (15) 

The density of states is thus D{to) ~ and the integrated density of the states N(u)) ~ in 
good agreement with what was found in the experiment and in the simulations. 

2.1. Comparision of methods 

In this section we briefly discuss and compare the two methods adopted in the present study 
to obtain the spectrum namely the 2x2 dynamical matrix K^y{k) and the spatial covariance 
matrix described in Eq. |6] To establish a correspondence we compute the full 2N x 2N matrix in 
k-space, 

Q^y(k,k')^{u;(k)uy(k')}. (16) 
using the same Fourier components of displacements as described earlier (Eq. |2), which in- 
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Figure 6: Normal mode spectrum: Cumulative Density of states and DOS for the present crystal as obtained from the 
covariance matrix 2pv(r,, Tj) in real space and its equivalent Q(k,k') in Fourier space. For a perfect crystal these two 
methods should give identical results. In the case of imperfections the modes \k) and \k') are no longer orthonormal 
leading to small differences in the spectra. 



eludes the non-orthogonal eorrelations as well. Let us first discuss the results obtained from the 
full diagonalization of both Q matrices (Eqs. l6land[T6Tl. Fig. |6]compares the spectrum - DOS 
and cumulative DOS for both the methods: they nearly coincide. This confirms the equivalence 
of the spatial matrix Q^vii, j) with Qfiv{k, k') and supports the fact that the data are not affected 
by the choice of the boundary conditions or imperfections of the crystal. Now, the diagonal ele- 
ments of Qijv(k, k') correspond to /r^! but the non-diagonal elements in {k, k') encode information 
about the heterogeneity and imperfections of the sample. Thus, with negligible imperfections in 
the crystal and a sufficiently long averaging time, the spectrum from all three methods should 
coincide with each other. 

2.2. Conclusions 

We have studied the dynamics of a hard sphere colloidal crystal at a volume fraction ~ 0.57, 
a volume fraction slightly above the melting transition, using confocal microscopy. The density 
of states and normal modes were obtained from measured particle displacements. Hard sphere 
systems are usually weakly connected and the interaction potential is strongly anharmonic; how- 
ever the present observations shows that the lowest frequency modes are extended plane waves- 
like as can be expected for a harmonic solid. In addition we have shown that the density of states 
can be understood using continuum elasticity theory. 

The effective exponent for the frequency-dependence of the density of states was measured 
in the low energy regime and is inconsistent with the expected Debye behavior in ~ w''"' 
for both d - 2 and d - i. We found that the data can be explained by a theory with an unusual 
energy dispersion relation in \k2\, which gives D{(jj) ~ aP. This expression agrees with both the 
experiments and the simulations. It is interesting to note that the same energy function was found 
in [16,1 where the spreading of a droplet was expressed as the effective dynamics of a contact line. 
Again we are in the presence of a physical system projected to lower dimensions. 

We appreciated discussions with Jorge Kurchan and Gerard Wegdam. The present project is 
supported by FOM. 
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